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NASTRAN.  The  solver  is  easily  modified  to  use  double  precision  complex  arithmetic 
on  computers  on  which  it  is  available.  With  somewhat  more  effort  It  could  be 
modified  to  solve  systems  of  equations  having  real-valued  coefficients. 
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Several  features  distinguish  this  linear  equation  solver  from  previous  solvers.  It 
automatically  determines  the  dimensions  of  the  blocks  of  coefficients  and  blocks  of 
right-hand  side  vectors  and  avoids  the  need  for  adding  extra  equations  by  providing  for 
blocks  of  more  than  one  size.  It  accepts  columns  of  the  coefficient  matrix  and  columns 
of  the  right-hand  side  vectors  from  a  sequential  file  and  returns  the  columns  of  solution 
vectors  on  a  sequential  file.  The  program  OCSOLVE  will  solve  with  one  call,  a  linear 
system  having  multiple  right-hand  sides.  It  will  solve  a  system  of  linear  equations  if 
a  specified  minimum  number  of  words  are  provided  for'  storing  the  blocks;  however,  the 
more  memory  provided,  and  hence  the  larger  the  blocks,  the  more  efficient  the  solution 
will  be. 
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ABSTRACT 

A  new  block  solver,  OCSOLVE,  for  large,  full,  unsymmetric 
systems  of  algebraic  equations  with  complex-valued  coefficients  has 
been  developed.  Although  OCSOLVE  was  developed  for  use  with  the 
finite  element  program  NASTRAN,  it  is  designed  to  be  easily  adapted 
for  other  applications.  This  new  solver  was  developed  because 
NASTRAN's  solver  was  not  designed  to  solve  full,  unsymmetric 
systems  efficiently;  it  reduced  the  time  required  to  solve  such  a  system 
of  500  equations  with  complex-valued  coefficients  to  about  5%  of  the 
time  required  by  the  equation  solver  in  NASTRAN.  The  solver  is 
easily  modified  to  use  double  precision  complex  arithmetic  on 
computers  on  which  it  is  available.  With  somewhat  more  effort  it 
could  be  modified  to  solve  systems  of  equations  having  real-valued 
coefficients. 

Several  features  distinguish  this  linear  equation  solver  from 
previous  solvers  It  automatically  determines  the  dimensions  of  the 
blocks  of  coefficients  and  blocks  of  right-hand  side  vectors  and  avoids 
the  need  for  adding  extra  equations  by  providing  for  blocks  of  more 
than  one  size  It  accepts  columns  of  the  coefficient  matrix  and 
columns  of  the  right-hand  side  vectors  from  a  sequential  file  and 
returns  the  columns  of  solution  vectors  on  a  sequential  file.  The 
program  OCSOLVE  will  solve  with  one  call,  a  linear  system  having 
multiple  right-hand  sides.  It  will  solve  a  system  of  lineal-  equations  if 
a  specified  minimum  number  of  words  are  provided  for  storing  the 
blocks;  however,  the  more  memory  provided,  and  hence  the  larger  the 
blocks,  the  more  efficient  the  solution  will  be. 


ADMINISTRATIVE  INFORMATION 

This  work  was  sponsored  by  the  Foundation  Acoustic  |)<-Mgn  Program.  Task  Area  S 1 255001  .  Element 
0356 9 N .  DTRC  Work  I  nit  l'Jll-s’01  The  Naval  S.-a  '■A-.tems  t'ommand  program  manager  was  R  ('Inn 
(NAVSEA  5553)  The  DTRC  program  manager  w  a-  (>  Ritter  Code  I 1  1 .  ami  the  te.  lineal  lin-uii, 
T  Tinley,  Code  1 7*20  1 


INTRODl  <  I  K  >N 


The  David  Taylor  Research  ('enter  |I>TR<  )  has  developed  a  capability  fur  computing  arotMie  pre-Min 
fields  outside  submerged  three-dimensional  eiasln  tl  n  >- 1 1 1 1  •  I  hi-  ■  apal  aiily  .  called  NAHM  \  1  -  l- 
llllple  tnelited  Using  the  finite  e  le  ni  e  n  1  -I'll'tillal  all  a  1  v  i  program  N\S!K\N  and  III'  'bides  a  lequip  III'  lit  I"! 
solving  a  large  system  of  linear  algebiaic  equation-  lie  ■  =  •  •  - 1 1 1  cM  iiiaUiV  bu  tin  v  -  t  <  •  1 1 1  is  lull  cninp|i\ 
it  lis  \  m  met  lie  and  so  large  that  the  iin'iiini  i  m  <  ■  b  ■  i  I '  ■  i  1<  uin  ■■  it  ■  \  ■  >  ■  I  tie  in  e  m  <  <1  \  n  ail  d  ,|e  >  n  the 
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Center’s  CDC  Cyber  176  and  Cray  X-MP.  Also  the  system  must  be  solved  for  many  rightrhand  sides.  To 


satisfy  the  requirements  of  NASHUA,  an  out-of-core  equation  solver  has  been  written  that  can  solve  this  type 
of  linear  system.  This  solver  is  designed  to  solve  a  system  of  equations  produced  by  another  program  (in  this 
case  NASTRAN)  that  writes  the  coefficients  and  rightrhand  side  vectors  on  a  sequential  file  and  requires  the 
solution  vectors  to  be  returned  on  a  sequential  file.  The  program  is  written  in  standard  FORTRAN  77  to 
facilitate  its  use  on  different  computers.  Although  the  program  is  written  for  single  precision  complex 
coefficients,  it  is  easily  modified  to  use  double  precision  complex  arithmetic  with  computers  for  which  this 
arithmetic  is  available.  Also,  but  with  somewhat  more  effort,  the  program  could  be  modified  to  solve  systems 
of  equations  having  real-valued  coefficients.  When  blocks  of  the  coefficient  matrix  and  of  the  right-hand  side 
vectors  are  in  central  memory,  they  are  stored  in  arrays  in  blank  COMMON.  For  single  precision  complex 
arithmetic,  the  minimum  memory  required  in  blank  COMMON  for  these  blocks  is  twice  the  length  of  one 
column  of  the  coefficient  matrix  plus  150  words.  Of  course,  the  larger  the  number  of  equations  and  the 
smaller  the  number  of  words  of  memory  provided  for  storing  the  blocks,  the  greater  the  cost  of  the  solution. 

These  systems  of  equations  can  be  solved  using  subroutines  provided  in  NASTRAN.  However,  since  it 
is  unusual  to  encounter  large,  full,  unsyminetric  systems  in  finite  element  structural  analyses,  the  procedures 
incorporated  in  NASTRAN  have  not  been  optimized  to  the  extent  that  have  the  more  commonly  used 
equation  -olvers  and  consequently  the  solution  is  very  slow.  For  this  application,  the  program  OCSOLVK 
reduced  the  time  to  solve  a  full  unsyni metric,  complex  system  of  order  500  to  about  5*7  of  the  time  required 
by  the  equation  solver  included  m  NASTRAN 
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Tliis  report  documents  the  out-ol-core  block  solver  program  OC’SOlAK  developed  by  the  Applied 
Mathematics  Division  at  l)TR('  to  solve  large  linear  systems  AX  B,  where  A  is  a  full,  u  iismii  met  i  ic  matrix 
of  complex  roe  tin  mnt.s,  li  is  a  matrix  of  one  or  more  right-hand  side  vectors,  and  X  is  the  matrix  of  unknown 
vrtois  I’his  program  solves  a  system  of  linear  equations  bv  partitioning  the  matrix  of  coc lli'  ient.s  and  the 
matrix  I  >1  right-hand  ide  Veilois  into  S||  bill  at  rices  called  blocks  The  blocks  ale  stored  III  direct  aeec.ss  |i|rs 
and  hiought  into  .  enti.nl  memory  a  lew  at  a  time  The  program  (X'SOIAK  will  solve  a  system  ol  any  older 
Willi  .111  V  ll  II  III  be  |  o|  I  I 1 1  I  -  1 1 .1 II  d  side  Vei  lei-  A  I  dor  k  solver  Is  most  e  1 1 1  r  |  r  n  t  il  the  lilorf-,  a|*  a-  large  a- 


possible  for  the  memory  available.  Therefore,  it  is  desirable  that  the  solver  automatically  determine  block 
sizes,  read  the  data  from  a  sequential  tile,  store  it  in  blocked  configuration,  and  after  solving  the  equations 
return  the  solution  vectors  in  a  sequential  file. 

Cantin'^  and  Hofmeister^  have  developed  block  solvers  for  large  systems  of  linear  equations.  C'antin's 
solver  is  designed  for  large,  symmetric,  banded  systems;  Hofmeister’s  solver  is  designed  for  large, 
unsymmetric,  full  systems.  Each  will  solve  a  system  of  equations  for  one  right-hand  side  vector  at  a  time. 
These  programs  require  that  the  matrices  have  already  been  stored  in  block  configuration  and  that  all  the 

blocks  are  of  equal  size.  To  satisfy  these  requirements,  the  user  must  organize  the  data  in  block  configuration 

and,  if  the  number  of  rows  of  blocks  does  not  divide  the  order  of  the  system,  must  add  equations  to  fill  out 
the  hist  row  and  column  of  blocks.  Care  is  needed  in  reorganizing  sequential  data  into  block  configuration  to 
avoid  excessive  calls  to  disk  storage,  for  although  one  access  to  disk  storage  is  fairly  quick,  many  calls  will  add 
significantly  to  the  time  required  for  solving  the  equations  If  all  blocks  must  be  of  equal  size  m  the  block 
elimination  process,  the  number  of  equations  added  to  the  system  will  be  no  greater  than  the  number  of  rows 
of  blocks  used  Thus  the  penalty  in  time  to  solve  the  equations  is  not  too  large  if  only  a  few  rows  o|  blocks 
are  used  However,  adding  the  extra  equations  require.',  another  'ti  p  m  the  prepiuaiion  of  the  system  for 
solution  and  requires  some  additional  time  m  the  'olutioii  <d  tin-  \stem  Providing  for  variably  sized  blocks 
eliminates  both  disadvantage'  (  an  tin  ml  llofm.  i-tei  both  us.-  m  a<  h  lie -depe  mb-iil  calls  to  access  disk 

storage  \  It  bough  it  i'  not  difficult  to  modify  tin  -•  i  ill  -  I  <  H I  P  \  N  77  m  w  pi  os  ides  standard  procedures  for 

reading  and  writing  to  direct  ace  e-  lib--  (  ic,ii\  ,  hoc  loped  .1  hick  d\ei  espec  ial!y  lor  systems  of 
equation'  arising  1 1 1  apple  ale  m  .  cl  i  h  •  I  u  n  d  u  v  i  n  n  c  i  ,  i  >  •  \  iiat  e .  u  in-  1 1 1  ■  i  III'  ■  ol \  e  r  ■  b  '  ign  ed  lor 

systems  With  bilge  blocks  of  /e|.i  *•  1 1 1 1  I  ill  ill'  •  ••til'  l  e  lit  lu.lUlV 

rile  block  I'll  III  ill  at  loll  pi".  •  "  II  '■'  I  III  ill'  pi  "Cl  INI  <  •<  '<  >1  \  I  I  '  •  -e  III  i.illy  tli.il  c  d  (  all  I  III  as  modi  lied 

b\  I  1 1  if  III  e  |-  ter  *  I  lie  program  I  III  pb  III '  II  t  III  c’  ill '  pi  o.  e  ■  fi.i.  I  ■  I :  •  I  ■  •  - 1  • ;  I ! '  I  I  •  male  il  c',1..'.  lei  II"'  ni'll  cctll'  I 

pi  ogi  am s  as  an  m le i  m ■■■  li  it<-  i ■  >ii t in •  I  In  i< df  •  w  m g  I .  . 1 1 . 1 1  •  a •  •  e ■  ■  w  t< ■  (  M  >|  \  | 

•  |l  solve-.  With  0,0  call  'Well!  id  equation  With  ■  1 1  I  i  I  •  1 1  I-  ll’lll  II  .III  •  I  id' 


•  It  obtains  the  order  of  the  system  and  the  number  of  right-hand  side  vectors  from  a  sequential  hie  and 
automatically  determines  a  suitable  blocking  of  the  matrix  of  coefficients  and  the  matrix  of  right-hand 
side  vectors. 

•  It  takes  these  matrices  from  the  sequential  file,  and  stores  them  in  the  block  configuration  required  by 
the  block  solver. 

•  It  reorganizes  the  solution  vectors  from  block  configuration  and  writes  them  on  a  sequential  file  after 
the  system  is  solved. 

•  It  does  not  require  the  blocks  to  have  equal  dimensions,  so  the  order  of  the  system  cf  equations  being 
solved  need  not  be  divisible  by  the  number  of  rows  of  blocks. 

The  program  OCSOLVE  is  designed  to  read  data  from  sequential  files  with  quite  general  format  and 
return  the  solution  in  a  sequential  file  with  an  equally  general  format.  To  do  this,  the  program  calls  live  very- 
short  subroutines  that  read  a  label,  heading,  or  matrix  column  from  the  data  file  or  write  a  heading  or  column 
to  the  solution  file.  Appendix  A  shows  an  example  of  a  set  of  these  routines  that  can  he  used  with  sequential 
files  with  a  simple,  fixed  format  and  another  set  that  can  be  used  with  sequential  files  with  the  NASTHAN 
OUTPUT2  or  INPUTT2  format.  If  the  solver  is  to  be  used  with  a  program  that  uses  and  produces  .sequential 
files  with  another  format,  the  user  can  easily  modify  these  routines  to  accommodate  that  format 

To  determine  the  size  of  the  submatrix  blocks,  the  program  OCSOLVE  uses  (lie  order  of  the  linear 
system  and  the  memory  available  for  the  blank  COMMON  arrays  The  order  of  the  linear  system  is  obtained 
from  the  datafile.  The  user  determines,  according  to  the  memory  available  in  the  particular  application,  tip- 
number  of  words  available  for  blank  COMMON  and  supplies  this  number  as  a  parameter  in  the  FORTRAN 
PARAMETER  statement  at  the  beginning  of  the  program.  The  user  also  provides  m  tlu>  statement  the  name 
of  the  file  that  contains  the  input  data  and  the  output  solution  vectors  as  well  as  parameters  ili.it  implement 
single  or  double  precision  arithmetic  and  indicate  the  type  of  computer  that  is  iped 


OUT-OF-CORE  STRATEGY 


In  the  out-of-core  elimination  strategy,  the  matrix  of  coefficients  and  the  matrix  of  rightrhand  side 
vectors  are  partitioned  into  blocks  of  submatrices,  producing  a  matrix-like  array  of  blocks.  The  block 
elimination  procedure  is  illustrated  here  for  a  system  of  four  equations. 

all  ai2  ai3  aH 

3 loi  3*22  2U>3  ^*24 
;l31  a32  a33  a34 

a4 1  a42  a43  a44 

The  coefficients  and  the  right-hand  side  terms  kj  are  known,  and  the  x-  are  to  be  determined  First  the 
solution  of  the  system  by  the  usual  Gaussian  elimination  will  be  described,  and  then  the  block  elimination 
procedure  will  be  shown  to  be  the  same  procedure  with  matrix  arithmetic  replacing  numeric  arithmetic. 

When  this  system  is  solved  by  Gaussian  elimination,^  elementary  row  operations  are  used  to  transform 
the  system  to  an  equivalent  system  in  which  all  entries  in  the  coefficient  matrix  below  the  main  diagonal  are 
zero.  For  this  example  the  equations  are 

bnx,  +  b|oX2  +  b13x3  b14x4  -  h, 

booX-j  ^  bo3-\ 3  b<2iX4  Ilo 

bsaXj  f  baix,4  h* 

t*44X 4  ~  ^4- 

The  last  equation  can  easily  be  solved  for  x4  b4.,1  h,  With  x,  known,  the  second  last  e<(tiatiou  can  be 

solved  for  x3  ■=-  b33l  (h3  -  b34x4).  Similarly,  each  of  the  remaining  two  equations  is  solved  In.  the  remaining 

unknown  values  x.j  and  x  1 ,  respectively.  The  elementary  row  operations  use  addition  and  multiple  ation  of  the 

numbers  a  and  k 
|J  i 

l,'or  the  block  elimination  procedure,  the  coefficient  matrix ,  the  unknown  vector,  and  the  right-hand 
side  vector  are  each  partitioned  to  form  two  matrix  equations  with  two  vector  unknowns  The  partition 
produces  a  total  of  eight  blocks  of  coefficients,  unknowns,  and  right-hand  side  valms 
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Elementary  row  operations,  now  using  multiplication  and  addition  of  matrices,  transform  the  system  into  the 
equivalent  system 
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The  lower  two  rows  form  a  linear  system  of  order  two  that  can  be  solved  by  the  usual  Gaussian  elimination 
for  x3  and  x4.  Substituting  these  values  into  the  upper  two  rows  produces  another  system  of  order  two  that 
can  be  solved  for  X)  and  x2.  In  this  example  of  a  block  solution,  only  systems  of  order  two  are  solved  and  the 
procedure  can  be  organized  so  that  at  most  four  blocks  are  required  in  memory  at  one  time.  When  not  in 
memory,  the  blocks  are  stored  in  direct  access  files.  Even  for  this  very  small  example,  an  in-core  solution 
would  require  that  all  eight  blocks  be  stored  in  memory  at  one  time.  The  blocks  are  assigned  to  records  in 
column  order  from  the  array  of  blocks.  The  matrix  entries  in  each  block  are  also  stored  in  column  order  in 
the  block’s  record  in  the  direct  access  file 

The  size  of  the  blocks  is  determined  by  the  order  of  the  linear  system,  the  number  of  right-hand  side 
vectors,  and  the  memory  available  If  the  order  of  the  matrix  of  eoellicients  is  divisible  by  the  number  of 
blocks,  the  coefficient  blocks  will  all  be  square  and  have  equal  dimensions,  otherwise  all  coellieient  blocks 
except,  those  in  the  Iasi  row  ami  column  will  be  square  and  have  equal  dimensions.  :uid  the  dimensions  of  the 
blocks  in  the  last  row  and  column  will  be  decreased  so  that  the  total  number  of  rows  and  columns  m  the  Mock 
configuration  equals  that  of  the  original  matrix  of  coefficients.  In  general  the  blocks  of  right-hand  side  vectors 
will  he  rectangular.  The  number  of  rows  in  a  right-hand  side  block  will  be  1  h same  as  tlm  numb. a  .,|  lows  m 


the  coellieient  blocks  m  the  same  block  row  The  numhei  o|  .  dunum  m  raefi  right-hand  id.  block  m 
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determined  by  the  size  of  the  problem  and  the  memory  available  for  storing  blocks. 

The  block  dimensions  are  computed  using  the  following  variable  definitions: 

NBRC  is  the  number  of  rows  and  columns  in  a  coefficient  block. 

NRHB  is  the  number  of  columns  in  a  rightrhand  side  vector  block. 

NRHV  is  the  number  of  right-hand  side  vectors. 

LBC  is  the  number  of  complex  numbers  available  in  blank  COMMON. 

At  times  there  will  be  in  memory  three  coefficient  blocks  of  dimension  NBRC  x  NBRC,  one  right-hand  side 
block  of  dimension  NBRC  x  NRHB,  one  column  of  working  storage,  and  one  column  of  pivots.  The  block  of 
right-hand  side  vectors  and  the  two  columns,  each  of  length  NBRC,  constitute  the  fourth  block.  The  first 
estimate  for  NBRC  follows  from  the  possibility  of  four  blocks  in  memory,  and  the  condition  that  the  fourth 
block  will  be  no  larger  than  the  first  three. 

4  x  NBRC2  =  LBC 

where  LBC  is  the  number  of  complex  numbers  that  can  be  contained  in  blank  COMMON.  The  number  of 
columns  in  a  rightrhand  side  block  is  determined  by  the  conditions  that  (1)  NRHB  is  at  most  NRIIV.  and  (2) 
the  fourth  block  in  the  preceding  estimate  includes  one  column  for  pivot  indices  and  one  column  for  working 


storage. 


NRHB  =  mini  NBRC-2,  NRHV 


The  program  reads  columns  of  the  coefficient  and  right-hand-side  matrices  into  the  fourth  block  while  storing 
the  system  in  block  configuration,  therefore,  this  block  must  he  at  least  large  enough  to  contain  one  column 
of  length  NMRC.  Thus,  the  number  of  columns  in  the  fourth  block  is 


NBK4  ■-  max(  NRHV  t- 2, 


NMRC  t  NBRC  1 
NBRC 


Then  the  equation 


3(  NBRC) “  t  BBKlxNBRC  LBC 


is  solved  to  determine  the  number  of  rows  and  columns  m  a  coefficient 


V  kV.v  .V '*Vr\V.%\VlV.V. V 


After  the  sizes  of  the  blocks  have  been  computed,  the  system  of  equations  is  stored  in  block 
configuration  on  direct  access  files.  Then  the  blocked  system  of  equations  is  solved  by  Gaussian  elimination 


The  block  elimination  procedure  follows  the  same  steps  as  conventional  Gaussian  elimination  for  linear 
3,6 

equations. 

ARRANGEMENT  OF  THE  PROGRAM  OCSOLYE 
'I'll e  program  OCSOLYE  consists  of  three  principal  subroutines,  three  secondary  subroutines,  and 
several  user-provided  subroutines.  The  principal  subroutines  are 

•  STORAB,  which  reads  the  matrix  of  coefficients  and  the  matrix  of  right-hand  side  vectors  and  stores 
them  in  block  configuration; 

•  BLKSLV,  which  performs  the  block  elimination  process;  and 

•  I  NPACK,  which  retrieves  the  solution  vectors  from  the  block  configuration  and  writes  them  on  a 
sequential  file . 


The  secondary  subroutines  are 


•  SOLVEC,  which  solves  a  linear  system  of  equations  with  complex  coefficients 


MI  LT,  which  multiplies  two  matrices;  and 


BLKSIZ,  which  computes  the  sizes  of  the  submatrix  blocks. 


The  remaining  user-provided  subroutines  read  and  write  the  data  files  (see  Appendix  A) 


I  SER-PROYIDED  PARAMETERS  AND  SUBROUTINES 


To  determine  the  size  of  the  blocks,  in  addition  to  the  sizes  of  the  matrices  of 


■tits  and  right- 


hand  side  vectors,  the  program  OOSOLVE  requires  information  on  the  memory  available  m  blank  <  (  >M  Ml  >N 
Blank  COMMON  contains  the  arrays  for  the  blocks  and  the  column  of  pivots,  and  thus  tin  -ize  of  blank 
COMMON  determines  how  large  the  (docks  may  be  lor -several  machine-dependent  p.itaiiutii  (  H  s(  )|  \| 
needs  the  tv  p e  o!  computer  being  Used  The  program  al-o  needs  the  logn  al  lile  name  a.- -l  au  d  to  tin  input 


f  ■r-  •/  s* 

V '-C  -V-N  -V-" 


\  V  v>  '  ,n-.'aV 
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and  output  tiles.  The  following  FORTRAN  statement  located  at  the  beginning  of  the  source  code  of  the 
program  OCSOLVE  provides  this  information: 


PARAMETER! MACHNR  =  m,LGTHBC=  n,D  ATFIL  - ‘name’.NPREC  i) 


The  integer  n  is  the  length  of  blank  COMMON  in  terms  of  real  words,  the  length  can  be  calculated  h\ 
subtracting  the  number  of  words  of  memory  occupied  bv  ;dl  coding  ;uid  variables  from  the  total  field 
length  (The  variable  LBC,  which  is  the  number  of  complex  numbers  in  blank  COMMON,  is 
computed  from  LCJTHBC  and  depends  oil  whether  single  or  double  precision  complex  numbers  are 
being  used.)  The  integer  in  is  a  parameter,  provided  by  the  user,  that  indicates  the  type  of  computer 
Used,  and  is  given  by 


COMPUTER 

CRAY 

CDC 

IBM 

VAX 

APOLLO 


MACH  NR 

m  1 
m  2 
m  -  3 
m  1 


The  character  string  ‘name’  is  the  logical  name  of  the  file  that  contains  the  matrix  of  coellicieni.s  ami 
the  right-hand  side  vectors  which  are  input  to  OCSOLVE,  and  the  solution  vectors  which  are  output 
from  OCSOLVE.  The  integer  i  is  set  bv  the  user  to  equal  either  1  or  2  when  implementing  single  or 
double  precision  arithmetic.  (When  changing  from  single  to  double  precision,  or  back,  a  lew  other 
changes  must  be  made;  see  Appendix  B.) 

I'lie  user-provided  subroutines  lead  from  and  write  to  (lie  data  and  solution  tile  I  lie  program 
OCSOLVE  reipi  ires  a  minimum  of  data  on  the  input  data  lib'  and  requires  the  data  to  he  in  a  p  re  sen  bed  order 
however,  the  file  may  contain  more  data  ami  the  specific  format  may  be  varied,  as  long  as  the  requirement.-  <d 
tins  section  are  -  iti-licd  The  input,  data  file,  111.1.  1  may  begin  with  a  label;  following  an v  label  mil  l  In  a 
beading  lfi.it  specifies  the  number  of  rows  and  columns  of  the  matrix  of  right-hand  side  lector-  l,o||owim.: 

tin-  fir-t  heading  are  the  records  containing  the  columns  <d  the  right-hand  side  vectors,  each  . . nl  ■  <  *  n  I  am  in ", 

at  mo-t  one  column  Next  n  a  second  heading  that  -.pecilies  the  number  n|  row-  and  column.  ■  •(  the  m  itiix 


of  coefficients  (equal  to  the  number  of  rows  of  the  right-hand  side  vectors)  followed  by  the  columns  of  the 
matrix  of  coefficients,  again  each  record  containing  at  most  one  column. 

The  data  file  is  read  by  three  user-provided  subroutines.  To  design  these  subroutines,  one  must  know 
the  format  of  the  data  file  provided  to  OCSOLVE  and  that  of  the  solution  file  to  be  returned.  Appendix  A 
gives  examples  of  the  user-provided  subroutines  designed  for  straight-forward  unformatted  FORTRAN  files 
and  for  unformatted  NASTRAN  files.  The  subroutines  must  be  designed  to  be  consistent  with  the  format  of 
the  data  files  and  also  must  satisfy  certain  specifications  to  be  consistent  with  OCSOLVE.  The  subroutines  and 
their  specifications  are 

•  RDLBL,  which  rewinds  the  data  file  and  then  reads  and  discards  any  information  between  the 
beginning  of  the  tile  and  the  heading.  If  there  is  no  information  before  the  heading,  this  subroutine 
must  at  least  rewind  FILE  4  The  subroutine  call  is 

GALL  RDLBL 

The  subroutine  RDLBL  must  leave  the  pointer  positioned  so  that  the  subroutine  RDIIDG  can  read  the 
heading 

•  RDHDG,  which  reads  the  heading  preceding  each  of  the  matrices  of  right-hand  side  veit<>r>  and 
coefficients.  On  the  first  call,  before  the  matrix  of  right-hand  side  vectors,  it  returns  NR.  the  number 
of  rows  of  a  right-hand  side  vector,  and  NO,  the  number  of  columns  of  rigbt-lmnd  -mb  ve.nus  on 
the  second  call,  before  the  matrix  of  eoellieient.s,  it  returns  NR  and  N<\  which  have  the  aim-  value 
and  are  equal  to  the  number  of  rows  and  columns  of  the  matrix  of  roellirjrnts.  The  suhtoutinc  cal!  i- 

CALL  RDM DG(  NR.NO) 

The  subroutine  RDHDG  must  leave  the  pointer  positioned  so  tli.it  the  -uhnmtme  R I ) (  <  >1  .  an  read  tin- 
lirst  column  of  the  matrix  that  follows 

•  RIK'OL.  which  reads  one  matrix  column  from  the  data  file  Tie-  subn.ulim  •  il  i 


GALL  RDGOLf  RIG'.NMRt  | 


where  RB2  is  the  array  into  which  the  column  is  read  and  NMRC  is  the  number  of  complex  numbers 
to  be  read.  After  each  of  these  calls  the  pointer  must  be  positioned  so  that  the  subroutine  RDCOL  can 
read  the  next  column  of  the  matrix.  After  the  last  column  of  rightrhand  side  vectors  has  been  read, 
the  pointer  must  be  positioned  so  that  the  subroutine  RDHDG  can  read  the  heading  before  the  matrix 
of  coefficients. 

After  the  system  of  equations  has  been  solved,  the  program  OCSOLVE  produces  a  sequential  file  with 
the  results  of  the  solution  of  the  linear  equations  This  file  contains  a  heading  that  may  contain  the  number 
of  rows  and  columns  of  solution  vectors,  followed  by  the  columns  of  solution  vectors.  Three  subroutines 
must  be  provided  that  are  designed  to  produce  the  output  file  in  the  format  desired  by  the  user  and  to 
conform  to  the  general  specifications  that  follow.  The  subroutine  WTHDG  is  called  once,  and  then  for  each 
solution  vector  the  subroutine  VVTCOL  is  called  once.  After  all  the  solution  vectors  have  been  written,  the 
subroutine  WTEND  is  called  once  to  allow  the  solution  file  to  be  closed  and  rewound. 

•  WTHDG,  which  writes  a  heading  on  the  output  file,  FILE  i  The  heading  may  include  NR  and  NO  as 
defined  above.  If  NR  and  NC  are  not  needed  by  the  program  receiving  the  solution  file,  the  heading 
need  not  be  written  on  the  solution  file,  b  a  subroutine  W'THDG  which  at  least  rewinds  FILE  I  must 
be  provided.  The  subroutine  call  is 

CALL  WTIIDGf  NR.NC) 

«  WTCOL,  which  writes  one  matrix  column  to  the  output  file  The  subroutine  call  is 
CALL  \VTCOL( AK.LSG2) 

where  Ah  is  the  array  from  which  the  column  is  written,  and  LSG2  is  the  number  of  complex  numbers 
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to  be  written. 

WTEND ,  which  writes  an  end-of-lile  mark  and  rewinds  the  output  file  The  subroutine  call  is 

CALL  WTEND 


1 1 


DISCUSSION 


The  block  solver  OCSOLVE  has  been  developed  to  solve  large,  full,  unsymmetric  systems  of  complex 
linear  equations.  It  has  been  used  as  an  intermediate  routine  in  an  analysis  in  which  it  received  data  from  a 
sequential  XASTRAN  output  file  and  returned  the  solution  in  a  sequential  NASTRAN  input  file.  Block 
solvers  have  been  developed  previously,  but  the  new  features  of  this  solver  are: 

•  It  determines  the  dimensions  of  the  blocks  of  coefficients  and  rightrhand  side  vectors. 

•  It  stores  the  coefficient  matrix  and  right-hand  side  vectors  in  blocked  form. 

•  It  avoids  the  need  for  extra  equations  by  providing  for  blocks  ol  different  sizes. 

•  It  solves  with  one  call  a  linear  system  of  any  order  having  any  number  of  right-hand  sides. 

•  It  stores  data  from  a  sequential  file  in  the  blocked  configuration. 

OC'SOI.VE  reduced  the  elapsed  time  in  a  dedicated  environment,  on  the  CDC  Cyber  176  computer  for  the 
solution  of  a  system  of  order  500  from  the  41  wall  clock  minutes  taken  by  the  equation  solver  in  NASTRAN, 
which  was  not  optimized  for  full,  complex,  unsymmetric  systems  of  equations,  to  approximately  'J  minutes. 

A  lme.ir  equation  solver  is  used  m  OCSOLVE  to  factor  liie  diagonal  blocks  and  reduce  the  off-diagonal 
blocks  This  solver  uses  partial  pivoting  to  reduce  roundoff  error  and  to  avoid  failure  if  a  zero  appears  on  the 
diagonal:  however,  the  search  for  the  largest  pivot  is  confined  to  the  block  in  memory,  and  it  is  possible  that, 
for  some  column,  all  pivot  candidates  in  that  block  would  be  inadequate  The  block  elimination  process  could 
be  modified  to  extend  the  pivot,  search  to  the  other  blocks  in  ihe  column  of  blocks,  but  the  modification 
would  be  expensive  to  develop  and  expensive  to  run  Also  the  equation  solver  used  does  not  equilibrate  the 
system  If  equilibration  is  desired,  the  subroutine  SOLVI'C  could  be  modified  or  replaced  with  a  solver  that 
does  equilibrate  the  system  It  is  expected  that,  lor  many  systems  ol  lineal  equations  arising  Irom  physical 
problems,  the  coefficients  would  be  of  similar  orders  of  magnitude,  and  there  would  be  Luge  entile-  on  , . i 
near  the  diagonal  Thus  there  would  be  no  difficulties  due  to  the  incomplete  patti.nl  pivoting  m  tin  bn  k  o| 
equilibration,  f  or  other  systems  of  linear  equations,  such  as  the  systems  with  blocks  ol  zeros  ifi.it  nine  when 
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APPENDIX  A.  EXAMPLES  OF  THE  USER-PROVIDED  SUBROUTINES 

TH is  appendix  provides  two  sets  of  subroutines  for  reading  the  data  file  that  contains  the  matrix  of  the 
right-hand  side  vectors  and  the  matrix  of  coefficients  for  the  system  of  linear  equations  to  he  solved  and  then 
writes  the  solution  matrix  on  the  same  file.  The  first  set  can  be  used  with  simple  unformatted  FOR 'IRAN 
files  and  the  second  with  unformatted  NASTRAN  OUTPUT2  and  INPUTT2  files  In  either  case  FILE  4  is 
opened  by  the  statement 

OPEN(4,FILE=DATFIL,FORM=’UNFORMATTED  ’) . 
where  the  character  string  DATFIL  has  been  set  by  a  FORTRAN  PARAMETER  statement  pee  \ppendix  B| 


Setl.  Unformatted  FORTRAN  Files 

Subroutines  from  this  set  can  be  used  to  read  unformatted  FORTRAN  data  tile-  which  contain  exartlv 
the  following  records: 

One  record  containing  in  two  integer  words  the  number  of  rows  and 
columns  of  the  right-hand  side  vectors, 

One  record  for  each  column  of  the  right-hand  side  vector, 

One  record  containing  in  two  integer  words  the  number  of  rows  and 
columns  of  the  coefficient  matrix,  and 

One  record  for  each  column  of  the  coefficient  matrix 

Subroutines  from  this  set  can  also  be  used  to  produce  an  unformatted  solution  file  which  contains  the  number 
of  rows  and  columns  of  the  solution  vectors  in  the  first  record  and  one  column  of  the  solution  matrix  m  each 
of  the  following  records.  The  subroutines  in  Set  (  are 

SUBROUTINE  RI)I, BE 
REWIND  l 
RETURN 
N 


■T*. 


SUBROUTINE  RDHDG(NR,NC) 

READ (4)  NR,NC 

RETURN 

END 


SUBROUTINE  RDCOL(  NMRC.RB) 
COMPLEX  RB(  *) 

READ(4)  (RB(I),  I  =  l.NMRC) 

RETURN 

END 


SUBROUTINE  WRTHDG(NR.NC) 

REWIND  4 

WRITE(4)  NR.NC 

RETURN 

END 


SUBROUTINE  VVTCOL(  AK.LSG2) 
COMPLEX  AK(  *) 

WRITE( 4)  ( AK(I),I  -  I,LSG2) 

RETURN 

END 


SUBROUTINE  WTEND 
END  FILE  4 
REWIND  1 
RETURN 


Set  2  NASTRAN  OUTPUT2  and  INPUTT2  Files 

Subroutines  from  thus  set  ran  he  used  to  lead  the  unformatted  NASTRAN  files  produced  by 
OUTPUT2  The  subroutines  RDLBL,  RDIIDG,  ami  RDCOL  are  provided  to  read  the  label,  heading,  and 
columns  of  a  NAS  TRAN  unformatted  OUTPUT2  file  Similarly  the  subroutines  W'TIIDG  and  WT(  Ol.  are 
provided  to  write  :ui  unformatted  NASTRAN  INl’f  TT2  tile*  The  subroutines  in  Set  2  are 

SUBROUTINE  RDLBL 
REWIND  l 
DO  100  I  l.S 
100  READ(  4) 

RETURN 

END 


SDBROt  TINi:  RDIIDG(NR.NU) 

COMMON  II DO  NI.N2.N;U\NTI  KNT2.I\NTUN<  L.NRW.IIORM  HAIM  N/W 
I  ID  I  NS  IT  LO .  NA  M  E(  2 ) 


1 1> 


THE  INPUT  AND  OUTPUT  FILES  ARE  STANDARD  NASTRAN  UTl  FILES,  WHERE 
THE  INPUT  FILE  HAS  A  LABEL  AND  THE  OUTPUT  FILE  DOES  NOT 
THE  SAME  FILE  (TAPE4=UTl)  IS  USED  FOR  BOTH  INPUT  AND  OUTPUT. 

ON  INPUT,  UTl  CONTAINS  THE  RHS  MATRIX  AND  THE  COEFFICIENT 
MATRIX  IN  THAT  ORDER  ON  OUTPUT,  UTl  CONTAINS  THE  SOLUTION 
MATRIX. 

THE  NASTRAN  DM AP  INSTRUCTIONS  TO  INTERFACE  WITH  ’OCSOLVE’  ARE 
OUTPUT2  RHS, COEF, .,//-!  $  TO  GENERATE  THE  MATRICES 
AND 

INPUTT2  /SOLN,,,,  $  TO  READ  THE  SOLUTION 

READ(4)  Nl 
READ(4)  NAME 
READ(  4)  KNTl 

READ) 4)  N2 

READ  (4)  KF ,  N  CL ,  N  R  W .  I  FOR  M ,  I  IT  PE .  N  Z  W ,  1 1)  ENS ,  I FL  G 
READ( 4)  KNT2 

READ(4)  N3 
READ(4)  NAME 
READ(4)  KNT3 

NR  —  NRW 
NC  NCI, 

RETURN 

END 

SUBROUTINE  WRTHDG 

COMMON  IIDG  N  U N2  N :i  K  Nil.  K  N  I  '_> .  K  NT.t  .NT  'L.NRW.I  IT  >R  M ,  I  TV  I  'I ; .  N  /.  W 
1  IDENS.II  Ld.NAMECJ) 

REWIND  t 

WRi  TE(  1)  NT 
WRITE!  1)  NAME 
WRITE!  I)  KNTl 

KI  RS  I  201 
WRITE!  4)  N2 

WRITE!  1)  KIRST.\(  I.  NRW  II  OHM  TIT  PI!  N/W  IDENS.II  IT. 

WRITE!  It  K N  12 


WRITI 

WRIT! 

WRITI 


1  I)  NT 
(  I)  NAME 
(  II  KNIT 


RETURN 

END 

SUBROUTINE  RIX  OI  (  NMR< '  RIO 


W.' 


COMPLEX  RB(*) 

READ(4)  IB 

READ(4)  (RB|I),  I  =  1,NMRC) 

READ(4)  KNT 

RETURN 

END 


SUBROUTINE  \VTCOL( AK.LSG2) 
COMPLEX  AK(  *) 

COMMON  /KOUNT/  KNT 

NMRC2  =  2*LSG2 

WRITE(  4)  NMRC2 

WRITE( 4)  (AK(I).I  =  1.LSG2) 

WRITE(  4)  KNT 

KNT  =  KNT- 1 

RETURN 

END 


SUBROUTINE  WTEND 
KNT  =  0 
WRITE(4)  KNT 
END  FILE  4 
REWIND  4 
RETURN 


/y 


APPENDIX  B.  USING  OCSOLVE 

The  program  OCSOLVE  must  be  prepared  by  setting  the  user-provided  parameters  and  including  the 
user-provided  subroutines  from  one  of  the  two  examples  in  Appendix  A,  or  by  providing  routines  that  satisfy 
the  requirements  of  Appendix  A.  Specifically,  the  following  steps  must  be  taken: 

1.  Provide  the  parameters  on  the  FORTRAN  PARAMETER  statement  located  at  the  beginning  of  the 


source  code  for  the  program  OCSOLVE. 

PARAMETER(MACHNR— m.LGTH  BC -•  n.l)  ATFI1, --‘name'.  NPREC  i) 

-  n  is  the  length  of  blank  COMMON  in  terms  of  real  words;  the  length  can  be 
calculated  by  subtracting  the  number  of  words  of  memory  occupied  by  til  1 
coding  and  variables  from  the  total  field  length 

-  m  is  a  parameter  that  indicates  the  type  of  computer  used,  and  is  given  by 

COMPUTER  MACHNR 


CRAY 

m 

=  i 

CDC 

in 

—  2 

IBM 

m 

=  3 

VAX 

m 

4 

APOLLO 

m 

---•  5 

-  'name’  is  a  character  string  that  is  is  the  logical  name  of  the  file  that 
contains  the  matrix  of  coefficients  and  the  right-hand  side  vectors  that 
are  input  to  OCSOLVE  and  contains  the  solution  vectors  that  are  output 
from  OCSOLVE. 

-  i  is  set  by  the  user  to  either  I  or  2  when  implementing  single  or  double 
precision  arithmetic  (  Also  see  additional  steps  gi\cn  below  for  changing 
to  double  precision  ) 


2  Provide  the  Sll  broil  lilies  R  1 )  I,  HI ,  HI  >  1 1  1  H  ;  .  and  1;  1  H  <  >1  to  read  the  data  tile  (see  a  ppell  dl\  A  ) 

,i  ProV  ide  the  subroutines  U  mix;,  W  I  (  '()(,  .  and  \\  II  .  \  It  to  u  I  lie  the  solution  (lie  (see  \  ppendix  A  | 

To  use  OCSOLVE  the  ll-e|  1 1 1  e  I  e  l  \  pro\  ides  t||e  lie.,  -ill  flat  a  all  d  input  tile-  and  then  exeeiltes 
OCSOLVE  The  -pel  |||C  Step-  me 


Provide  a  data  file  with  the  logical  name  'perilled  b\  the  pi  e.  .duo!  I '  \  I .’  \  \  l  I  II  1 - 1  il<  no  ni 


Execute  (  X  )|.\  1 


I  'I 


3.  Get  the  solution  from  the  file  that  has  the  same  logical  name  as  the  input  file. 


If  the  program  OCSOLVE  encounters  a  singular  diagonal  block,  it  will  stop  with  the  message: 

SINGULARITY  IN  KNN  IN  POSITION  n 

where  n  is  the  number  of  the  equation  for  which  the  solver  could  not  obtain  a  nonzero  pivot. 

To  change  from  single  to  double  precision  complex  arithmetic,  the  parameter  i  on  the  PARAMETER 
statement  must  be  set  to  2,  and  all  TYPE  COMPLEX  declaration  statements  must  be  changed  to  TYPE 
COMPLEXM6  statements. 
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